Performance of Quantum Heat Engines Enhanced by Adiabatic Deformation of Trapping Potential

We present a quantum Otto engine model alternatively driven by a hot and a cold heat reservoir and consisting of two isochoric and two adiabatic strokes, where the adiabatic expansion or compression is realized by adiabatically changing the shape of the potential. Here, we show that such an adiabatic deformation may alter operation mode and enhance machine performance by increasing output work and efficiency, even with the advantage of decreasing work fluctuations. If the heat engine in the sudden limit operates under maximal power by optimizing the control parameter, the efficiency shows certain universal behavior, η*=ηC/2+ηC2/8+O(ηC3), where ηC=1−βhr/βcr is the Carnot efficiency, with βhr(βcr) being the inverse temperature of the hot (cold) reservoir. However, such efficiency under maximal power can be produced by our machine model in the regimes where the machine without adiabatic deformation can only operate as a heater or a refrigerator.

To describe the machine performance, there are usually two benchmark parameters: [2,8,17,27,32]: the thermodynamic efficiency η = w / q h , where w is the average work output per cycle and q h is the average heat released from the hot reservoir, and the power P = w /τ cyc , with the cycle period τ cyc . Ideally, both these two quantities should have large values for excellent performance, but there is always a power-efficiency trade-off dilemma [43][44][45][46][47][48][49]. An important issue is, hence, that of optimizing the heat engines by ensuring their efficiency under maximal power [1,2,7,11,20,50,51].
Discreteness of energy levels, due to quantization, may significantly improve the performance of a quasi-static quantum Otto cycle [19,31,32,[52][53][54] when an inhomogeneous shift of energy levels occurs along an isentropic, adiabatic stroke [54,55]. However, the question as to how such a shift (due to adiabatic deformation of potential) affects a quantum heat engine in the finite-time cycle period, as hinted at in [54], has not previously been answered. Moreover, the random transitions between discrete energy levels are responsible for quantum fluctuations, which dominate at low enough temperatures. A question naturally arises: what is the influence of such adiabatic deformation of potential, related to discrete energy spectra, on the relative power fluctuations that measure engine stability?
As we demonstrate, machine efficiency can be improved via controlling the shape of the potential, without sacrificing of machine stability.
In this paper, we study a quantum version of the Otto engine model, which consists of a single particle confined in two different potentials and which works between two heat reservoirs of constant inverse temperatures β r h and β r c (>β r h ). We analyze the engine performance by determining the efficiency and power with respect to the times spent on the two isochoric strokes and adiabatic deformation parameters. Assuming that only the two lowest energy levels are populated, we show that the adiabatic shape deformation of the potential operates as a heat engine in regions where its counterpart, without the deformation, works as a heater or a refrigerator. We find that the adiabatic deformation enhances performance, as well as stability, through appropriately selecting the forms of the two trapping potentials. We highlight that, in the sudden limit where the total time spent on the two adiabatic strokes is negligible, the efficiency at maximum power of our model shows universal behavior: Yet, such optimized efficiency can be obtained in the regions where the machine, in the absence of adiabatic deformation, cannot operate as a heat engine.

A Single Particle in a Power-Law Trap
We consider a single particle with mass m confined in a one-dimensional power-law trapping potential V(x) along x direction with This simple class of traps covers, for instance, harmonic (θ = 3/2), spherical-quadrupole (θ = 3), and infinite potential (θ = 0) traps. The Hamiltonian system,Ĥ =p 2 /2m + V(x), with momentum operatorp, can be written in terms of a single particle energy spectrum ε n , whereâ † n (â n ) is the creation (annihilation) operator, with single-particle quantum number n. Thus,â † nân is the particle number operator with quantum number n. Throughout the paper we seth ≡ 1 for simplicity. The energy spectrum can be written as ε n = n|Ĥ|n = ωn σ . Here, ω is the energy gap between the ground state and the first excited state (which we call the energy gap, for simplicity, in what follows), and σ (>0) is called the potential exponent [56], which is determined by the parameter θ in Equation (1) (For the onedimensional trapped system in the x direction, the time-independent Schrödinger equation may be written as [p 2 /2m + V(x)]Ψ n (x) = ε n Ψ n (x), where ε n are the energy eigenvalues and Ψ n are eigenfunctions. For example, for a single particle confined in a box trap which reads V(x) = 0 for 0 ≤ x ≤ L and V(x) = ∞ [due to θ = 0 with θ defined in Equation (1)] otherwise, the energy spectrum is obtained as ε n = n 2 π 2 /(2mL 2 ). When the single particle is confined in a harmonic trap with θ = 2/3, the potential becomes V(x) = ωx 2 /2, leading to ε n = (1/2 + n)ω. The energy eigenvalues ε n = ωn σ are determined by the trapping potential V(x), and the shape of trapping potential associated with θ can be captured well by the so-called trap exponent σ.) and is, thus, dependent on the shape of the external potential. For example, for a one-dimensional harmonic trap σ = 1 and ε n = nω , where ω is the trap frequency, and for a one-dimensional infinite deep potential (also called a one-dimensional box trap) with length L, σ = 2 and ε n = n 2 ω, where ω ≡ π 2 /(2mL 2 ).
The expressions for creation and annihilation operators (â † n andâ n ) in Equation (2) depend on the trapping potential under consideration. A typical example is that a single particle is confined in a one-dimensional infinite potential well, which is given by V(x) = 0 for 0 ≤ x ≤ L and V(x) = ∞ otherwise. The state wave function of a trapped particle reads Ψ n (s) = |n = √ 2/L sin(ns) with s = πx/L when 0 ≤ s ≤ π. The creation and annihilation operators for the system should satisfyâ † n |n ∼ |n + 1 andâ n |n ∼ |n − 1 [57]. In view of the fact that dψ n (s)/ds = √ 2/L ncos(ns), we define the creation and annihilation operators asâ n = Here,N is the number operator defined byN|n = √ n|n and its inverseN −1 |n = n −1 |n . By using these definitions, we obtainâ n |n = √ n|n − 1 and a † n |n = √ n + 1|n + 1 . We then find that the commutator is [â n ,â † n ] = 1 , and the system energy becomes Ĥ = ∑ n n 2 /(2mL 2 ) â † nân , where â † nân corresponds to the occupation probability at state n for the single-particle system.
The state of the system at thermal equilibrium with a heat bath of inverse temperature β can be described by the canonical formρ = ∑ n p n |n n| = Z −1 exp(−βĤ), where p n = e −βε n /Z is the probability of finding the system in state |n , with the partition function Z = Tr(e −βĤ ). The system entropy reads S = −Tr(ρ lnρ), whereρ =ρ(βω, σ), and, thus, the entropy takes the form of S = S(βω, σ). For the gas in a given trap, the entropy S merely depends on the parameter βω: S = S(βω), and in an adiabatic process βω = constant. However, an adiabatic deformation of trap, by changing σ, leads to change in the parameter 'βω' [54,55] to keep entropy S constant. That is, a quantum adiabatic process where the entropy is kept constant can be realized via changing the shape of the trapping potential.

General Expressions of Efficiency and Power for Quantum Otto Engines with Deformation of Trapping Potential
In contrast to conventional quantum heat engines, where the working substance is confined in a given form of trap, the quantum engine under consideration works, based on two different forms of one-dimensional trapping potentials V(x), by adopting two different values of θ in Equation (1). The quantum Otto engine, sketched in Figure 1, consists of four consecutive strokes, as outlined in the following: (i) Hot isochoric stroke A → B. The single particle is confined in a one-dimensional trap along the x direction with θ = θ h , and the trapped system is weakly coupled to a hot reservoir of constant inverse temperature β r h in time duration τ h . Since the external field V(x) is frozen, the energy gap is kept constant at ω = ω h ; (ii) Adiabatic expansion B → C. The von Neumann entropy of the system is constant along the adiabatic stroke in which the system evolution is unitary. While the system is isolated from the heat reservoir in time τ hc , and the form of the potential gradually changes from the trap V(x) ∼ x 3/θ h to the trap V(x) ∼ x 3/θ c by tuning θ; (iii) Cold isochoric stroke C → D. Both the trap configuration and the energy gap are kept fixed, namely, ω = ω c and θ = θ c . Within a time interval of τ c , the trapped system is weakly coupled to a cold reservoir with constant inverse temperature β r c (> β r h ); (iv) Adiabatic compression D → A. The system is again isolated from the heat reservoir in time duration τ ch , and the trap configuration changes gradually from the trap V(x) ∼ x 3/θ c to the trap V(x) ∼ x 3/θ h . During the hot or cold isochoric strokes, the system would relax to the thermal state at the ending instant B(D) of the hot (cold) isochore, if τ h (τ c ) is long enough . The times allocated to the four strokes set the total cycle period, For the Otto cycle, the work is produced only in the two adiabatic branches, with heat produced alongside the isochoric processes. Initially, the time is assumed to be t = 0. The Hamiltonian system changes fromĤ(τ h ) toĤ(τ h + τ hc ) along the adiabatic expansion B → C, and it goes back toĤ(0) fromĤ(τ cyc − τ ch ) after the adiabatic compression D → A. The Hamiltonian system is kept constant along each isochoric stroke, namely, . The stochastic work done by the system, per cycle, is, thus, the total work output along the two adiabatic trajectories [32,58] The stochastic work for the engine cycle is then given by where we used During the adiabatic stroke, the level populations do not change, p n,B = p n,C and p m,A = p m,D , and the probability density of the stochastic work w can then be determined according to where δ(·) is the Dirac's δ function. The average work output per cycle, w = wp(w)dw, can be obtained as Here, and hereafter, we use the subscripts A, B, C, and D (in Figure 1) to indicate the physical quantity at times t = 0, τ h , τ h + τ hc , and t = τ cyc − τ c , respectively. We define the dimensionless energy g as g ≡ Tr(ρĤ)/(ω), which, according to Equation (2), can be expressed as As emphasized,â † nân indicates the occupation number operator of a given state n, and, thus, â † nân is the average occupation number at state n. While the trapping potentials are V(x) ∼ x 3/θ h and V(x) ∼ x 3/θ c in the hot and cold isochoric strokes, respectively, we use the trap exponents σ h and σ c (rather than θ h and θ c ) to characterize the forms of the trapping potential in what follows. The engine model works with a single particle confined in the trap. It consists of two isochoric processes A → B and C → D, where the system is weakly coupled to the hot and cold heat reservoirs of constant inverse temperatures β r h and β r c (> β r h ), respectively, and two adiabatic processes B → C and D → A, where the shape of the trap is adiabatically deformed by changing θ in Equation (1) from θ = θ h to θ = θ c , or vice versa. In each cycle, the average work output ( w ) is the difference between the heat absorbed from the hot bath ( q h ) and the heat released to the cold reservoir ( q c ); that is, To describe the degree of the shape deformation of the trapping potential, we introduce the deformation parameters for adiabatic compression and expansion which are defined by respectively. Except in the special case when the shape of the potential is not changed along the engine cycle with σ c = σ h and ξ hc = ξ ch = 1, these parameters, ξ hc and ξ ch , depend on the trap exponents σ c and σ h , and, thus, they capture all information about the adiabatic deformation of trapping potential. We also note that, in the presence of adiabatic deformation, the deformation parameters, ξ ch and ξ hc , would be affected by the times, τ c and τ h , since the so-called system temperatures, at instants A, B, C, and D in Figure 1, are dependent on these times τ c and τ h . Using Equations (5)- (7), we find that the average work takes the form of The work fluctuations can be determined according to where Based on the two-time measurement approach, the probability density function of the stochastic heat q h along the hot isochoric stroke, where no work is produced, can be determined by the conditional probability to arrive at where p m,A is the probability that the system is initially in state m at time t = 0, and p τ h m→n is the probability of the system collapsing into another state n after a time period τ h . Here p For each cycle, heat is transferred only in the isochore, while work is produced only along the adiabatic process. The heat absorbed from the hot bath is given by Due to energy conservation, the heat discharged to the cold reservoir along the cold isochoric stroke can be directly calculated according to q c = q h − w (see also Figure 1). In order to evaluate the average values of heat and work in a finite time cycle, we should analyze the system dynamics along two isochoric strokes to derive average work and heat. We use Γ c (Γ h ) to denote the thermal conductivity between the system and cold (hot) heat reservoir and introduce x = e −Γ h τ h and y = e −Γ c τ c . We show that these quantities, (8) and (11), can be expressed as a function of x and y (see Appendix A for details), and where we used G = (1−x)(1−y) (1−ξ hc ξ ch xy) . The heat quantity released into the cold bath can be directly calculated by q c = w − q h , due to the conservation of energy. In the absence of adiabatic deformation, these average values, (12) and (13), reduce to w = (ω h − ω c ) g eq h − g eq c G, and q h = ω h g eq h − g eq c G. In such a case, we present these formulae in a broader context by considering a power-law trap in which g may not be the mean population if σ c,h = 1. We reproduce the result obtained from the harmonic trap, where σ c = σ h = 1, and, thus, g denotes the average population.
The thermodynamic efficiency, η = w / q h , then follows as which In the case when the shape of the potential is adiabatically changed, an inhomogeneous shift of energy levels is created, resulting in thermodynamic efficiency (14) that depends on the shapes of the potentials along two isochoric strokes, excepting the case when the two potentials are identical to each other, which would result in efficiency reducing to that of the cycles without adiabatic shape deformation, η = 1 − ω c /ω h .

Performance and Stability of a Two-Level Machine
The efficiency may be enhanced by adiabatically changing the form of the potential. To better understand the influence induced by adiabatic deformation on the performance of thermal machine, we investigate how the adiabatic deformation affects the efficiency and the power. In this section, we consider, as an example, the Otto engine working in the low-temperature limit, by assuming that only the two lowest energy levels are appreciably populated. We show in Appendix B that, for the two-level engine where the system Hamiltonian (2) simplifies toĤ = ω(â † 1â 1 + 2 σâ † 2â 2 ) and the system energy becomes Ĥ = ω ∑ n=1,2 n σ â † nân , the deformation parameters ξ ch and ξ hc defined by Equation (7) take the forms of where g eq c = e −β r c ωc +γ c e −γc β r c ωc Equation (15) into Equations (12) and (14), it follows that the average work (12) and thermodynamic efficiency (14) of the two-level machine in finite time are given by and Note that when γ c < γ h , the efficiency (17) is larger than the efficiency without adiabatic deformation (η = 1 − ω c /ω h ). For the two-level engine, the work fluctuation δ 2 w (9) can be analytically obtained as where g B and g D are When adiabatic deformation is absent, the work fluctuations turn out to be The system reaches thermal equilibrium at the end of the hot or cold isochore when the process is in the quasi-static limit. In this case, where x → 0, y → 0, and G = (1−x)(1−y) 1−xy → 1, the work (16) and work fluctuations (18) turn out to be The power output and power fluctuations are then determined according to P = w /τ cyc and δ 2 P = δ 2 w /τ 2 cyc . In Figure 2a we plot the normalized efficiency η/η C at the quasi-static limit as a function of the ratio r (with r ≡ √ ω h /ω c ) in the presence of adiabatic shape deformation, comparing the corresponding result for the Otto engine without deformation of trap. In the absence of adiabatic deformation of trap (γ c = γ h ), the three different conditions of the compression ratio r correspond to the three modes of the machine: (1) for r ≤ 1, the machine operates as a heater, (2) for 1 < r ≤ r C ≡ β r c /β r h , it works as a heat engine, and (3) for r > r C , it becomes a refrigerator. However, when adiabatically changing the shape of the trapping potential, the machine can operate as a heat engine even in boundaries (1) and (3).  For complete thermalization along each isochore, Figure 3a displays the average work (19) and the work fluctuations (20) as a function of the compression ratio r, respectively, for γ c = γ h = 2, 2γ c = γ h = 4, and γ c = 2γ h = 4. While the efficiency is improved by increasing r for given γ c and γ h , the average work w first increases, and then decreases as the ratio r increases. The behavior of curves for work fluctuations δ 2 w as a function of r is dependent on γ c and γ h . It can be observed from Figure 3a that, while for γ c ≤ γ h the curve of work fluctuations, as a function of r, is linear, it becomes parabolic when γ c > γ h . Both the work fluctuations δ 2 w and average work w for γ c > γ h are much smaller than those obtained from the case when γ c ≤ γ h . Figure 3a also shows that the regime of positive work ( w > 0) is sensitively dependent on the values of the parameters γ c and γ h , and the presence of adiabatic deformation changes the positive work condition for the quantum engine. The coefficient of variation for power f P = √ δ 2 P /P, equivalent to the square root of the relative work fluctuations, δ 2 w / w , is also called the relative power fluctuation. This coefficient measures the dispersion of the probability distribution and, thus, can describe the machine stability [17]. Comparing the efficiency and relative power fluctuations of the engine with γ c = 2 (γ c = 4) to each other, Figure 3b shows that optimization of the quantum heat engine can be realized by selecting the appropriate form of trapping potential during the hot isochoric stroke. For example, the engine with γ c = γ h = 4 works at efficiency η = 0.76 and relative power fluctuation f P = 25.6, but the model with γ c = 4 and γ h = 3.6 operates under η = 0.73 and f P = 12.7. That is, when γ c = 4, the relative power fluctuation f p for γ h = 3.6 is halved as compared to its value in the absence of adiabatic deformation, while efficiency only slightly decreased. Another example to consider γ c = 2 is comparing γ h = 2 with γ h = 1.84. In contrast to the former case, where η = 0.76 and f P = 13.77, in the latter case η = 0.72 and f P = 6.9, showing again that the relative power fluctuations can be significantly decreased with a particularly small decrease in efficiency. More importantly, by suitably choosing the shapes of the trapping potential, we may even design an engine model that runs more stably and effectively, as also shown in Figure 3b. A typical example is that of the engine of γ c = 4 and γ h = 3.35, producing efficiency η = 0.7 with f P = 9.15, but the model with γ c = 2 and γ h = 1.9 runs at η = 0.74 and f P = 8.53. By comparison, the latter model shows better overall performance than the former, since it runs more stably by decreasing relative power fluctuations, even with higher efficiency η. Hence, adiabatically changing the form of the trapping potential may even contribute to a decrease in the relative power fluctuations with an increase in efficiency, when compared with an engine without adiabatic deformation. So, for quantum engines, having selected suitable machine parameters, adiabatic deformation along an engine cycle may be an effective optimization approach to significantly improve engine stability.
Both Figures 2 and 3 show that the average work, work fluctuations, coefficient of variation for power, and even positive work condition, are strongly affected by change in the values of γ c and γ h . To further see clearly how the adiabatic deformation quantitatively affects the performance and fluctuations for the engine, we plotted average work, work fluctuations, and coefficient of variation for power as a function of Carnot efficiency in Figure 4a-c, respectively, where the value of γ c was kept fixed (γ c = 2), and the value of γ h was slightly changed (γ h = 2, 1.8, 2.2). Figure 4a demonstrates that, in the positive work region, adiabatic deformation could increase the average work in the certain regime of η C , though it may decrease the work when η C was relatively large. The fluctuations of the heat engine, including the work fluctuations and the relative work fluctuations (coefficient of variation for power), were always significantly decreased by the adiabatic deformation [as per Figure 4b,c]. These figures show that the shape change of the trap may enhance the average work, unless the difference between the two bath reservoirs is particularly large, while it always enhances the machine stability captured by the fluctuations. As a specific example, at η C = 0.85, both the work and efficiency [see Equation (17)] could be enlarged, but the fluctuations of work and relative power decreased. In the present case, the adiabatic deformation was, thus, an essential ingredient in improving the engine performance and stability in a certain regime.

The Engine under Maximal Power Output
Since the power output, P = w /τ cyc , would vanish if the cycle was quasi-static and the cycle period approached infinity, the engine should operate, practically, in finite time to produce finite power output. In this section we consider the efficiency and power statistics for the two-level machine under maximum power by optimizing power with respect to external degrees, and assuming that the two adiabatic strokes are realized in the sudden limit [59,60]. It is not difficult to verify that the efficiency at maximum power η * can be determined by using the method shown in Appendix B to analytically obtain [11] , which shows the same universality with the CA efficiency [1,3,48] Figure 5 we plotted the analytical efficiency of maximum power η * anal as a function of the η C , comparing the exact numerical result for different values of γ c,h and the CA efficiency η CA . These curves of the optimal efficiency for different γ c,h , together with the analytical expression of η * anal , collapse into a single line, and they are in nice agreement with the CA efficiency η CA . It was, therefore, shown that the efficiency at maximum power, agreeing well with η CA , is independent of the shapes of the two trapping potentials. As emphasized, the heat engine under consideration can proceed at such efficiency in regimes where the machine, without adiabatic deformation, may only operate as a heater or a refrigerator. Figure 5. Plots of analytical expression η * anal and exact numerical calculations for efficiency at maximum power, and plot of the CA efficiency η CA . We use η (γ c ,γ h ) to denote the exact values of optimal efficiency for given γ c and γ h . The inverse temperature of hot bath was β r h = 2.

Conclusions
As a result of energy quantization, a quantum adiabatic process can be realized by changing the shape of the trapping potential. Such a shape deformation causes the classical limit, where the principle of the equipartition of energy holds, to vanish and is, therefore, of purely quantum origin. Here, we investigated the performance and fluctuations in quantum Otto engines in the presence of adiabatic deformation of one power-trap potential. We started using stochastic thermodynamics and the quantum master equation to determine heat and work statistics, and then presented general expressions for time-dependent efficiency and work [cf. Equations (12) and (14)], in which adiabatic deformation parameters of the two adiabatic strokes are involved.
We proposed an exact analytical description for the performance and fluctuations in these quantum engines at the low-temperature limit where only the lowest two energy levels are occupied. We showed that quantum heat engines, with adiabatic deformation, can run in the extended regimes where their counterparts, without adiabatic deformation, operate as heaters or refrigerators. Examining the efficiency and coefficient of variation of power, we found that an appropriate selection of two trapping potentials enables engines to be built that are capable of performing more stably and efficiently. We also showed that, even for a given trap in an isochore, the relative power fluctuations in our engines are significantly smaller than those of engines in the absence of adiabatic deformation, with higher efficiency than that of engines without adiabatic deformation. By tuning the energy gap between the ground and the first excited state, we found that the efficiency at maximum power is independent of shape deformation and shares the same universality with the CA efficiency. This optimized efficiency, however, can be realized in regions where engines experiencing no change in the shape of the potentials cannot operate as heat engines.
Our approach can be directly used to describe an ensemble of many non-interacting particles (with particle number N) confined in a d−dimensional power-law trap. In such a case, the Hamiltonian system becomesĤ = ∑ n ε nâ † nân ,where ε n = n|Ĥ|n , with n = n 1 , · · · , n d and N = ∑ nâ † nân . Following the approach adopted in this paper, we can reproduce the same results and arrive at the same conclusions. Our results significantly add to the study of quantum heat engines in finite time by taking advantage of adiabatic deformation of trapping potential, facilitating the design of efficient and stable quantum heat engines.  with their Hermitian conjugatesV † α , and k α are phenomenological positive coefficients. SubstitutingĤ =X into Equation (A1), the first law of thermodynamics is obtained as, where w is stochastic work and q is the stochastic heat. Thus, dw/dt = ∂Ĥ/∂t and dq/dt = L D (Ĥ) are the identified instantaneous average power and average heat flux, respectively. The Kraus operatorsV † andV are chosen asV = ∑ nân andV † = ∑ nâ † n . By insertingX =Ĥ [with the Hamiltonian system (2)] into Equation (A2) and using [â n ,â † m ] = δ mn , we then obtain L D (Ĥ) = ∑ n ε n k u Here, the phenomenological coefficients k u n and k d n can be parameterized as k u n =ñ th Γ and k d n = (ñ th + 1)Γ, whereñ th = 1/(e βε n − 1) is the excitation number corresponding to the n-th energy level and Γ, as a constant, is associated with the systembath interaction strength. The detailed balance, k u n /k d n = e −βε n , is satisfied and the system can be driven to reach a thermal state after an infinitely long time. It then follows that, in the isochore with static Hamiltonian, Equation (A1) becomes where Ĥ eq ≡ ∑ n ε n k u n /(k d n − k u n ) = ∑ n ε nñth is the asymptotic value of Ĥ at thermal equilibrium, and Γ denotes the heat conductance between the system and reservoir. Since dĤ/dt = d Ĥ /dt, where Ĥ = ωg with g defined by Equation (6) in the main text, Equation (A3) can be simplified to where we used g eq = Ĥ eq /(ω). Now we are in a position to discuss the evolution of the system during the system-bath interaction interval. For the finite-time process, the system, assumed to be initially at time t A = 0, evolves from the initial instant A to the final state B. On this branch, heat is absorbed from the hot bath during a period τ h while no work is done. The system would relax to the thermal state after an infinitely long time, and then g B τ h →∞ = g eq h = g(β r h ω h , σ h ). From Equation (A4), we obtain where Γ h is the heat conductivity between the working substance and the hot reservoir. For the cold isochore C → D, the system is in contact with the cold reservoir at inverse temperature β r c for a time of τ c . Based on an analogy with the hot isochore A → B, the dimensionless system energy g D , as a function g C , is obtained, where g eq c = g(β r c ω c , σ c ), and Γ c represents the heat conductivity between the working substance and the cold reservoir. Having these formulae [Equations (7), (A5) and (A6)], the relationship between g(τ h ) [g(τ cyc − τ ch )] and its asymptotic value g eq h (g eq c ) is easily obtained: (1 − ξ hc ξ ch y)x 1 − ξ hc ξ ch xy + g eq c ξ ch (1 − y)x 1 − ξ hc ξ ch xy , g D = g eq c − g eq c (1 − ξ hc ξ ch x)y 1 − ξ hc ξ ch xy + g eq h ξ hc (1 − x)y 1 − ξ hc ξ ch xy , where we used x = e −Γ h τ h and y = e −Γ c τ c . As expected, x, y → 0, g B τ h →∞ → g eq h , g D τ c →∞ → g eq c when τ c,h → ∞. Inserting Equation (A7) into Equations (8) and (11) in the main text leads to the time-dependent expressions of work and heat injection [cf. Equations (12) and (13)].

Appendix B. Work (Power) Statistics for a Two-Level System
For an ideal gas confined in a one-dimensional power-law potential, the single-particle energy spectrum takes the form ε n = ωn σ , where n = 1, 2, · · · , and σ is a positive index depending on the form of the trapping potential. Several special case examples include the following [56]: (i) σ = 2 for an infinite potential well; (ii) σ = 1 for a harmonic potential; (iii) σ = 4/3 for a quartic potential. When only the first two levels are populated for onedimensional potential, the probabilities of these two levels at the end of the cold isochore are given by where Z D = e −β D ω c + e −γ c β D ω c with γ c = 2 σ c , and these two-level probabilities at the end of the hot isochore become where Z B = e −β B ω h + e −γ h β B ω h with γ h = 2 σ h . For example, when the harmonic trap in the cold isochore to the one-dimensional infinite potential well, then σ c = 1 and σ h = 2, but if the potential is one-dimensional box (harmonic) in the cold (hot) isochore, then σ c = 2 and σ h = 1. In such a case, the dimensionless energies g [cf. Equation (6)] at the four instants [in Figure 2a] can be written as where we used χ c = e −(γ c −1)β D ω c and χ h = e −(γ h −1)β B ω h . Here, γ c = γ h , if, and only if, σ c = σ h . Note that there is a relation: At the times τ c,h → ∞, the system reaches thermal equilibrium at the end of either the hot or the cold isochore, indicating that β B → β r h and β D → β r c . We, therefore, obtain g Combining Equation (A11) with Equations (A5), (A6) and Equation (7) , we then obtain Equation (15) in the main text. Inserting Equation (15) into Equations (12) and (14), we obtain the expressions of average work and efficiency [Equations (16) and (17) in the main text], as well as the power output